Quaternions

So gimbals, 3 accumulated axial rotations, do not really work very well for orienting an object. Gimbals can be locked, and it is very unintuitive to control them. How do we fix these problems?

Part of the problem is that we are trying to store an orientation as a series of 3 accumulated axial rotations. Orientations are orientations, not rotations. And orientations are certainly not a series of rotations. So we need to treat the orientation of the ship as an orientation, as a specific quantity.

The first thought towards this end would be to keep the orientation as a matrix. When the time comes to modify the orientation, we simply apply a transformation to this matrix, storing the result as the new current orientation.

This means that every yaw, pitch, and roll applied to the current orientation will be relative to that current orientation. Which is precisely what we need. If the user applies a positive yaw, you want that yaw to rotate them relative to where they are current pointing, not relative to some fixed coordinate system.

There are a few downsides to this approach. First, a 4x4 matrix is rather larger than 3 floating-point angles. But a much more difficult issue is that successive floating-point math can lead to errors. If you keep accumulating successive transformations of an object, once every 1/30th of a second for a period of several minutes or hours, these floating-point errors start accumulating. Eventually, the orientation stops being a pure rotation and starts incorporating scale and skewing characteristics.

The solution here is to re-orthonormalize the matrix after applying each transform. A coordinate system (which a matrix defines) is said to be orthonormal if the basis vectors are of unit length (no scale) and each axis is perpendicular to all of the others.

Unfortunately, re-orthonormalizing a matrix is not a simple operation. You could try to normalize each of the axis vectors with typical vector normalization, but that would not ensure that the matrix was orthonormal. It would remove scaling, but the axes would not be guaranteed to be perpendicular.

Orthonormalization is certainly possible. But there are better solutions. Such as using something called a quaternion.

A quaternion is (for the purposes of this conversation) a 4-dimensional vector that is treated in a special way. Any pure orientation change from one coordinate system to another can be represented by a rotation about some axis by some angle. A quaternion is a way of encoding this angle/axis rotation:

Equation 8.1. Angle/Axis to Quaternion

Axis = x y z Angle = θ Quaternion = x * sin θ 2 y * sin θ 2 z * sin θ 2 cos θ 2

Assuming the axis itself is a unit vector, this will produce a unit quaternion. That is, a quaternion with a length of 1.

Quaternions can be considered to be two parts: a vector part and a scalar part. The vector part are the first three components, when displayed in the order above. The scalar part is the last part.

Quaternion Math

Quaternions are equivalent to orientation matrices. You can compose two orientation quaternions using a special operation called quaternion multiplication. Given the quaternions a and b, the product of them is:

Equation 8.2. Quaternion Multiplication

a * b = a w * b x + a x * b w + a y * b z - a z * b y a w * b y + a y * b w + a z * b x - a x * b z a w * b z + a z * b w + a x * b y - a y * b x a w * b w - a x * b x - a y * b y - a z * b z

If the two quaternions being multiplied represent orientations, then the product of them is a composite orientation. This works like matrix multiplication, except only for orientations. Like matrix multiplication, quaternion multiplication is associative ( (a*b) * c = a * (b*c) ), but not commutative ( a*b != b*a ).

The main difference between matrices and quaternions that matters for our needs is that it is easy to keep a quaternion normalized. Simply perform a vector normalization on it after every few multiplications. This enables us to add numerous small rotations together without numerical precision problems showing up.

There is one more thing we need to be able to do: convert a quaternion into a rotation matrix. While we could convert a unit quaternion back into angle/axis rotations, it's much preferable to do it directly:

Equation 8.3. Quaternion to Matrix

1 - 2 q y q y - 2 q z q z 2 q x q y - 2 q w q z 2 q x q z + 2 q w q y 0 2 q x q y + 2 q w q z 1 - 2 q x q x - 2 q z q z 2 q y q z - 2 q w q x 0 2 q x q z - 2 q w q y 2 q y q z + 2 q w q x 1 - 2 q x q x - 2 q y q y 0 0 0 0 1

This does look suspiciously similar to the formula for generating a matrix from an angle/axis rotation.

Composition Type

So our goal is to compose successive rotations into a final orientation. When we want to increase the pitch, for example, we will take the current orientation and multiply into it a quaternion that represents a pitch rotation of a few degrees. The result becomes the new orientation.

But which side do we do the multiplication on? Quaternion multiplication is not commutative, so this will have an affect on the output. Well, it works exactly like matrix math.

Our positions (p) are in model space. We are transforming them into world space. The current transformation matrix is represented by the orientation O. Thus, to transform points, we use O*p

Now, we want to adjust the orientation O by applying some small pitch change. Well, the pitch of the model is defined by model space. Therefore, the pitch change (R) is a transformation that takes coordinates in model space and transforms them to the pitch space. So our total transformation is O*R*p ; the new orientation is O*R .

Yaw Pitch Roll

We implement this in the Quaternion YPR tutorial. This tutorial does not show gimbals, but the same controls exist for yaw, pitch, and roll transformations. Here, pressing the SpaceBar will switch between right-multiplying the YPR values to the current orientation and left-multiplying them. Post-multiplication will apply the YPR transforms from world-space.

Figure 8.3. Quaternion YPR Project

Quaternion YPR Project

The rendering code is pretty straightforward.

Example 8.2. Quaternion YPR Display

void display()
{
    glClearColor(0.0f, 0.0f, 0.0f, 0.0f);
    glClearDepth(1.0f);
    glClear(GL_COLOR_BUFFER_BIT | GL_DEPTH_BUFFER_BIT);
    
    glutil::MatrixStack currMatrix;
    currMatrix.Translate(glm::vec3(0.0f, 0.0f, -200.0f));
    currMatrix.ApplyMatrix(glm::mat4_cast(g_orientation));
    
    glUseProgram(theProgram);
    currMatrix.Scale(3.0, 3.0, 3.0);
    currMatrix.RotateX(-90);
    //Set the base color for this object.
    glUniform4f(baseColorUnif, 1.0, 1.0, 1.0, 1.0);
    glUniformMatrix4fv(modelToCameraMatrixUnif, 1, GL_FALSE, glm::value_ptr(currMatrix.Top()));
    
    g_pShip->Render("tint");
    
    glUseProgram(0);
    
    glutSwapBuffers();
}

Though GLSL does not have quaternion types or quaternion arithmetic, the GLM math library provides both. The g_orientation variable is of the type glm::fquat, which is a floating-point quaternion. The glm::mat4_cast function converts a quaternion into a 4x4 rotation matrix. This stands in place of the series of 3 rotations used in the last tutorial.

In response to keypresses, g_orientation is modified, applying a transform to it. This is done with the OffsetOrientation function.

Example 8.3. OffsetOrientation Function

void OffsetOrientation(const glm::vec3 &_axis, float fAngDeg)
{
    float fAngRad = Framework::DegToRad(fAngDeg);
    
    glm::vec3 axis = glm::normalize(_axis);
    
    axis = axis * sinf(fAngRad / 2.0f);
    float scalar = cosf(fAngRad / 2.0f);
    
    glm::fquat offset(scalar, axis.x, axis.y, axis.z);
    
    if(g_bRightMultiply)
        g_orientation = g_orientation * offset;
    else
        g_orientation = offset * g_orientation;
    
    g_orientation = glm::normalize(g_orientation);
}

This generates the offset quaternion from an angle and axis. Since the axis is normalized, there is no need to normalize the resulting offset quaternion. Then the offset is multiplied into the orientation, and the result is normalized.

In particular, pay attention to the difference between right multiplication and left multiplication. When you right-multiply, the offset orientation is in model space. When you left-multiply, the offset is in world space. Both of these can be useful for different purposes.

Fork me on GitHub